gusucode.com > 7波段参数均衡器源码程序 > 7波段参数均衡器源码程序/code/Parametric EQ/lowpassEQ.m

    function [Hd]=lowpassEQ(fc,g,N)

fs=44100;

% parameter limits 
if fc < 20
    fc=20;
end
if fc > 20000
    fc=20000;
end
if g < -12
    g=-12;
end
if g > 6;
    g=6;
end

% coefficient calculation
if N==1
    g=10^(g/20);
    a1=-0.5*tan(pi*(fc/fs-0.25));
    b0=g*(0.25-0.5*a1);
    b1=b0;
    b2=0;
    a2=0;
    % SOS Form
    b=2*[b0 b1 b2];
    a=[1 -2*a1 -2*a2];
    [sos,g]=tf2sos(b,a);
    Hd=dfilt.df1sos(sos,g);
elseif N==2
    g=10^(g/20);
    t0=tan(pi*fc/fs);
    t1=1+sqrt(2)*t0+t0^2;
    a1=(1-t0^2)/t1;
    a2=0.5*(sqrt(2)*t0-1-t0^2)/t1;
    b0=g*t0^2/(2*t1);
    b1=2*b0;
    b2=b0;
    % SOS Form
    b=2*[b0 b1 b2];
    a=[1 -2*a1 -2*a2];
    [sos,g]=tf2sos(b,a);
    Hd=dfilt.df1sos(sos,g);
elseif N==3
    g=10^(g/20);
    t0=tan(pi*fc/fs);
    t1=1+t0+t0^2;
    a01=0.5*(1-t0)/(1+t0);
    a11=(1-t0^2)/t1;
    a12=(t0-t0^2-1)/(2*t1);
    b00=g/2*1/(t0+1);
    b01=b00;
    b10=t0^2/(2*t1);
    b11=2*b10;
    b12=b10;
    a02=0;
    b02=0;
    % SOS Form
    b0=2*[b00 b01 b02];
    a0=[1 -2*a01 -2*a02];
    [sos0,g0]=tf2sos(b0,a0);
    Hd0=dfilt.df1sos(sos0,g0);
    b1=2*[b10 b11 b12];
    a1=[1 -2*a11 -2*a12];
    [sos1,g1]=tf2sos(b1,a1);
    Hd1=dfilt.df1sos(sos1,g1);
    Hd=dfilt.cascade(Hd0,Hd1);
elseif N==4
    g=10^(g/20);
    t0=tan(pi*fc/fs);
    t1=1+2*cos(1/8*pi)*t0+t0^2;
    t2=1+2*cos(3/8*pi)*t0+t0^2;    
    a01=(1-t0^2)/t1;
    a02=(2*cos(1/8*pi)*t0-1-t0^2)/(2*t1);
    a11=(1-t0^2)/t2;
    a12=(2*cos(3/8*pi)*t0-1-t0^2)/(2*t2);
    b00=0.5*t0^2/t1;
    b01=2*b00;
    b02=b00;
    b10=g*t0^2/(2*t2);
    b11=2*b10;
    b12=b10;
    % SOS Form
    b0=2*[b00 b01 b02];
    a0=[1 -2*a01 -2*a02];
    [sos0,g0]=tf2sos(b0,a0);
    Hd0=dfilt.df1sos(sos0,g0);
    b1=2*[b10 b11 b12];
    a1=[1 -2*a11 -2*a12];
    [sos1,g1]=tf2sos(b1,a1);
    Hd1=dfilt.df1sos(sos1,g1);
    Hd=dfilt.cascade(Hd0,Hd1);
end


if nargout==0
    freq=0:fs/4095:22050;
    h=zeros(4096,1);
    h(1)=1;
    h=filter(Hd,h);
    H=20*log10(abs(fft(h,4096)));
    fplot(freq,H(1:end/2));
end